Aging in attraction-driven colloidal glasses 
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Aging in an attraction-driven colloidal glass is studied by computer simulations. The system is 
equilibrated without attraction and instantaneously "quenched" , at constant colloid volume fraction, 
to one of two states beyond the glass transition; one is close to the transition, and the other one deep 
in the glass. The evolution of structural properties shows that bonds form in the system, increasing 
the local density, creating density deficits (holes) elsewhere. This process slows down with the time 
elapsed since the quench. As a consequence of bond formation, there is a slowing down of the 
dynamics, as measured by the mean squared displacement and the density, bond, and environment 
correlation functions. The density correlations can be time-rescaled to collapse their long time 
(structural) decay. The time scale for structural relaxation shows for both quenches a super-linear 
dependence on waiting time; it grows faster than the bond lifetime, showing the collective origin of 
the transition. At long waiting times and high attraction strength, we observe completely arrested 
dynamics for more than three decades in time, although individual bonds are not permanent on 
this time scale. The localization length decreases as the state moves deeper in the glass; the non- 
ergodicity parameter oscillates in phase with the structure factor. Our main results are obtained 
for systems with a barrier in the pair potential that inhibits phase separation. However, when this 
barrier is removed for the case of a deep quench, we find changes in the static structure but almost 
none in the dynamics. Hence our results for the aging behavior remain relevant to experiments in 
which the glass transition competes with phase separation. 

PACS numbers: 82.70.Dd, 64.70.Pf, 82.70.Gg 



I. INTRODUCTION 

When a liquid is quenched below its glass transition 
temperature, its dynamics slows down enormously and 
the structural relaxation time becomes bigger than the 
experimental time: an amorphous solid is formed. In 
structural glasses, in contrast to spin glasses where theo- 
retical approaches have been developed 0, no theory is 
yet available to describe the slowing down of states be- 
yond the glass transition, i.e., aging. Computer simula- 
tions have provided, nevertheless, important information 
about general properties of aging 0- 

In colloidal systems with short range attractions, vit- 
rification can be repulsion driven (induced by increasing 
the pressure) or attraction driven (induced by decreas- 
ing temperature at intermediate densities) [H, |j, IB] • The 
former is caused by the steric hindrance of hard cores at 
high densities (as is the glass transition in most atomic or 
molecular systems). The latter is more specific to colloids 
with short range attractions. It is induced by increasing 
the attraction strength, and caused by bonding between 
the particles. At high attraction strength and high den- 
sity, where both glass transitions merge, a higher order 
singularity appears, causing a logarithmic decay of the 
density correlation function |a B • Mode coupling 
theory (MCT) anticipated this scenario [Tol ITTl IT ^ . and 
correctly describes the liquid states close to both glass 
transitions, where the structural slowing down sets in 
0, H, [13, EJ, Il3 ■ However, it does not provide informa- 
tion about aging within the arrested state. 



Using molecular simulations, Kob et al studied the ag- 
ing of a structural (repulsive) glass of Lennard- Jones par- 
ticles . Quenching from high temperatures, the slowing 
down of the dynamics was studied by means of a two-time 
density correlation function, <$> q (t,t w ), where the earlier 
of the two times, t w , represents a 'waiting time' after the 
quench. They found a scaling of the structural relaxation 
(similar to the time-temperature scaling in fluid states) 
with the relaxation time following a power law with t w . 
More recent studies within the nonergodic region have fo- 
cused on dynamical hetero gen eities and their importance 
to the glass transition [H, 12 0] • 

Recent simulations of the attraction driven glass have 
considered systems where structural arrest interplays 
with liquid-gas separation, producing states of arrested 
phase separation 0, [53, Ej • The observed phenomena 
from simulations resemble the experimental observations 
|22| , but the dynamics of the system was not analyzed in 
terms of glass aging. In this paper we isolate the effects 
of the attractive glass transition from those of phase sep- 
aration, and study the resulting aging dynamics without 
this additional complication. One paper has addressed 
the attractive glass by studying a system at high concen- 
trations, above the onset of the single-phase liquid [23| . 
Starting from the fluid pocket between the attractive and 
repulsive glasses, the repulsive glass could be reached by 
increasing the temperature, and the attractive glass by 
decreasing it. However, the states studied are close to 
the higher order singularity mentioned above; this affects 
the dynamics near the glass transition 0, and possibly 
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also its aging within the non-ergodic regions. A different 
work by the same group |24j studied the transition from 
repulsive glass to attractive glass at high density, and it 
was concluded that the bonds were not strong enough 
stabilize the attractive glass and the glass to glass tran- 
sition was not observed. Introducing an artificial barrier 
to stabilize these bonds, arrested dynamics due to the at- 
tractive glass was finally obtained and the glass to glass 
transition observed p4| . 

In this work, we study a system which we have used 
previously to analyze the fluid close to the attraction 
driven glass transition HEHH. The liquid- gas transi- 
tion has been inhibited by introducing a long-range re- 
pulsive barrier in the interaction potential @ ■ We study 
the aging of the system in two different states, on an iso- 
chore where we have previously determined the transition 
point. On performing an instantaneous quench from the 
system without attraction, the structure and dynamics 
of the system are monitored. As bonds between particles 
form, the local density increases, causing holes to open 
up elsewhere in the system. Accordingly, structural re- 
laxation takes longer and longer as aging proceeds. 

We show below that the attractive glass is broadly sim- 
ilar to the repulsive one (and thus to atomic glasses) in 
many fundamental properties, such as the time-scaling 
of correlation functions and the evolution of the localiza- 
tion length. An important difference is also noticed be- 
tween the two states studied here: the glass close to the 
transition is "softer" than that deeper in the non-ergodic 
region, where the dynamics indeed arrests. The bond 
correlation function shows that the bonds are reversible 
even for the strongest attraction studied. The bond and 
density correlation functions evolve differently, in con- 
tradiction with recent results found in a different sys- 
tem at h igh temperatures and by a different simulation 
method |2q . Finally, we turn off the repulsive barrier 
and study the competition between the attractive glass 
transition and liquid-gas separation, to elucidate the role 
of the barrier in the dynamics. This investigation is mo- 
tivated in part by recent observations of equilibrium clus- 
ters in protein solutions at lower packing fractions where 
a charge induced repulsive barrier is considered essential 
[27j • The results show that the dynamics of the system 
is similar with and without the barrier. This suggests 
that the aging dynamics is primarily controlled by bond 
formation in a state of arrested phase separation, rather 
than slow evolution of the long-wavelength density fluc- 
tuations that are present in such a state. A companion 
study of dynamic heterogeneities in the attractive glass 
will be presented elsewhere |28|. 



II. SIMULATION DETAILS 

The simulated system consisted of 1000 quasi-hard par- 
ticles interacting by a short range attraction, with New- 
tonian dynamics. The system is polydisperse to avoid 
crystallization: particle sizes are distributed according 



to a flat distribution of half-width S = 0.1a, where a is 
the mean radius. The core-core repulsion is given by 



V, (t) = k B T 
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where 2ai2 = (ai + a 2 ), with a\ and a 2 the radii of the 
interacting particles. The short range attraction mim- 
ics the interaction between colloidal particles induced by 
non-adsorbing polymers in a colloid polymer mixture. 
For a polydisperse colloidal system, the interaction is 



VAo(r) = -k B Tcj> p 
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for 2ai2 < r < 2(ai2+£) and for larger distances. Here, 
r}i = aj/^; f] = (r/ 1 + 772 )/2, and 4> p is the volume fraction 
of the polymer. The range of the interaction, £, is the 
polymer size, and its strength is proportional to 4> p , the 
concentration of ideal polymers. This potential has been 
slightly corrected near r = 2a\ 2 to ensure (2^| that the 
total interaction potential has its minimum at 2ai2- 

At high polymer fractions, or attraction strengths, this 
system would undergo liquid-gas separation into two fluid 
phases. The interplay of this demixing with the attrac- 
tive glass transition is complex 22], and would compli- 
cate the interpretation of our results. To avoid it, we add 
a repulsive long-range barrier to the interaction potential: 



Vbar{r) =k B T ■ 



, [ ^ 4 -2f^-^V + l) (3) 
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for tq < r < r 2 and zero otherwise, with t\ = (r 2 +ro)/2. 
The limits of the barrier were set to tq = 2(ai2 + £), 
and r 2 = 4:a, and its height is only lkgT (equal to the 
attraction strength at 4> p = 0.0625). The barrier raises 
the energy of a dense phase, so that liquid-gas separation 
is suppressed. The resulting total interaction, V to t = 
V sc + Vao + Vbar, is analytical everywhere and allows 
straightforward integration of the equations of motion. 
A short range attraction complemented by long range 
repulsive may induce microphase separation [301 131| . We 
have performed some simulations without the repulsive 
barrier to test its effects on the dynamics of the system, 
as shown in Sec. III.C. 

Lengths are measured in units of the average radius, a; 
time in units of \J 4a 2 /3v 2 , where the thermal velocity v 
is set to -\/4/3. Equations of motion are integrated using 
the velocity- Verlet algorithm, in the canonical ensemble 
(constant NTV), to mimic the colloidal dynamics. Every 
n t time steps, the velocity of the particles is re-scaled to 
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assure constant temperature. The time step was set to 
0.0025 and n t = 100. The range of the attraction, 2£, 
is set to 2£ = 0.2. The density of colloids is reported as 

volume fraction, cf> c = ^1 + (§)~^ n c with n c the 

colloid number density; the attraction strength is fixed 
by the polymer volume fraction, <p p . (Note that, with 
£ = 0.1, the attractive minimum is, for average sized 
particles, at r = 2a, with V m i n = — 16fcsT0 p .) 

In our simulations the system is first equilibrated at 
4> p = and 4> c = 0.40 (infinite temperature) with the 
repulsive barrier, and then instantaneously "quenched" , 
at constant colloid density, to the desired attraction 
strength. The time elapsed since the quench is called the 
waiting time, t w . Both the static and dynamic quantities 
depend on t w ; i.e., structural properties become waiting- 
time dependent, and dynamical ones are two-time func- 
tions, depending on t w (the first of the two times) and 
t' = t — t w , with t the second time and t' the interval. 
For instance, the density correlation function reads, 
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FIG. 1: Evolution of the mean number of neighbors per 
particle (continuous lines, left scale) and energy per particle 
(broken lines, right scale) for both quenches: 4> p = 0.50 (black 
lines) and 4> p — 0.80 (red lines) .The blue horizontal line in the 
upper right corner marks the mean number of neighbors for 
the equilibrium state <f> p = 0.42. 



® q (t',t w ) = — ^exp{iq(rj(t' + t UJ ) ~r k (t w ))} (4) 

where j and k run over all N particles in the system. The 
mean squared displacement, the self part of the density 
correlation function, or bond correlation functions, can 
be defined similarly. 

The system is quenched to two states inside the non- 
ergodic region, <f) p = 0.50 and <j> p = 0.80 (both with col- 
loid density <fr c — 0.40). For each state, 50 independent 
simulations were run, to improve the ensemble average. 
In each simulation different waiting times were studied 
and the total duration of the simulation was t = 3x 10 4 . 
In the last part of this work (Sec. III.C), the repulsive 
barrier is switched off to study the competition between 
phase separation and the attractive glass transition in 
the state <p p — 0.80 (50 independent simulations were 
also performed for this case). 



III. RESULTS AND DISCUSSION 

Upon increasing the attraction strength at constant 
colloid density (<f> c = 0.40) the attractive glass transition 
is found at cj> p = 0.4265, as estimated from the analysis 
of the density correlation functions, self-diffusion coeffi- 
cients and viscosity according to MCT 0, [2f|, [3^] . 

As stated above, we study the aging dynamics for 
two states quenched to within the non-ergodic region, 
4> p = 0.50 and <p p — 0.80. The former is just within, 
and the latter well within, the non-ergodic region. (Note 
that the minimum in the potential is V m i n = —&ksT and 
V m i n = — 12.8/cbT for <f) p — 0.50 and (j) p = 0.80, respec- 
tively, for the average sized particles.) The structure and 
dynamics is tackled for both quenches and compared with 
an ergodic fluid close to the transition, namely 4> p = 0.42. 



Some theoretical predictions are also tested. The effect 
of the repulsive barrier is studied in the final subsection. 



A. Structure 

The instantaneous quench from <j) p = provokes the 
formation of bonds between particles. Since both states 
are above the percolation line in density, a bonded cluster 
extending throughout the system is formed at zero time 
(using the interaction range, 2£ = 0.2a, to define the 
length of a bond). 

In Fig. ^the evolution of the mean number of neigh- 
bors (bonds) per particle, and of the energy per particle, 
are presented for both quenches. As the bonds between 
particles form, the energy per particle decreases. This 
process slows down as the waiting time increases, but a 
steady state is not reached in the simulation. Interest- 
ingly, more bonds are formed in the quench to the lower 
4> p , although the energy per bond is smaller in this case. 
However, the equilibrium system at <f> p = 0.42 reaches a 
larger number of bonds per particle (horizontal blue line 
in the upper right corner of the figure), at least in the 
time studied in the simulation. Note that, in compari- 
son with repulsion driven glasses (with broadly similar 
potential minima) [2|, the energy evolves much more in 
the attractive case. This can be expected given the rela- 
tively low density [<f> c = 0.40) of our system. Recall that, 
at <j) p = no bonds are present in the system; but as 
soon as the quench is performed a highly bonded state 
exists. This state evolves with time, and can alter con- 
siderably since particles are not strongly caged by their 
hard-core repulsions at this modest <p c . This leaves much 
more room for energy drift in attractive glasses than in 
the repulsion-driven case. 

The intricate structure of the network can induce a 
heterogeneous distribution of the bonds per particle. In 
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FIG. 2: Distribution of bonds per particle for (j>p = 0.50 



and 



0.80 at different waiting times, as labeled. The 



discontinuous line shows the distribution in the equilibrium 
system at cf> p — 0.42. 



Fig. [21 the evolution of the distribution of bonds per par- 
ticle is presented for both states at different waiting times 
and compared with the distribution for the equilibrium 
system at 4> p — 0.42. All of the distributions (includ- 
ing that of the system in equilibrium) are wide, due to 
the heterogeneity of the system, but single-peaked. As 
aging proceeds, the distribution moves to higher num- 
bers of neighbors per particle causing an increase in the 
mean number, but never reaches the distribution at equi- 
librium. Thus, the difference in the mean number of 
neighbors shown in Fig. ^ is caused by a displacement 
of the distribution towards larger number of neighbors, 
and not by the development of stronger heterogeneities 
in the structure of the system. 

The structure factor of the system is presented in Fig. 
[3J for both quenches, at different waiting times. Increased 
local order is noticed in the main peak in S(q); this 
slightly increases with waiting time, and the oscillations 
at high-g become stronger. But the most prominent fea- 
ture in the evolution of S(q) is the appearance of a low-q 
peak, which grows and moves to smaller wavevectors with 
waiting time (see insets). Since the system has bond per- 
colation from the outset, this peak is associated not with 
cluster formation, but with voids appearing and growing 
(thus the peak moves to smaller q) as the bonding state 
evolves. Accordingly, the structure factor for the homo- 
geneous system at 4> p = shows no peak in this region, 




FIG. 3: Structure factor for different waiting times (as la- 
beled), for 4> p = 0.50 (upper panel) and <j> p — 0.80 (lower 
panel). S(q) for 4> p — and <j> p = 0.42 are also presented 
(broken and continuous thin green lines, respectively). The 
insets show the low-g region for every state magnified. 



and the neighbor peak is displaced to larger wavevectors 
in both quenches (particles come closer to each other). 

In agreement with Fig. ^ the fluid system at 4> p = 
0.42 shows a higher neighbor peak and stronger oscil- 
lations, since there are more bonds in the equilibrium 
state than in the quenches (the difference being larger 
for <f) p = 0.80). Also, the low-g peak is highest in the 
equilibrium state (on the time scales studied in the simu- 
lation), and higher for lower <fi p in the non-ergodic region. 

The overall trend is consistent with slow restructur- 
ing of the system to reach a lower energy; for stronger 
bonds the system is more hindered on its way towards 
the energy minimum (even though this is deeper in prin- 
ciple). This hindrance is characteristic of aging systems, 
although in attractive glasses it is caused by bonding, as 
opposed to core-core repulsions in repulsive glasses. 

The low wavevector peak in S(q) resembles a low an- 
gle peak found in experiments on colloidal gelation and 
aggregation [^U 02 03] • The origin of the experimen- 
tal peak has not fully been determined; since gelation 
inhibits liquid-gas separation, the peak is probably the 
signature of frustrated demixing (2J, • It also obeys 
the factorization rule typical of spinodal decomposition 
0ESHil In our case, the low q peak is not caused by 
an arrested phase separation, but it is an intrinsic charac- 
teristic of the system induced by the long range repulsive 
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FIG. 4: Mean squared displacement as a function of t' = t—t w 
for both quenches at different waiting times: from left to right 
t w = 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192 and 
16384. The thin green lines show the mean squared displace- 
ments for 4> p = (broken line) and <f> p = 0.42 (continuous 
line) . 

barrier in the interaction. It is a signature of the trend 
to microphase separation. Therefore, it will stop grow- 
ing and moving to smaller q since the system does not 
demix. Its behavior is thus similar to the experimental 
one but not of the same origin. We show below that the 
dynamics of the system is independent of the origin of 
this peak in our simulations. 

B. Dynamics 

The stiffness of the system with increasing attraction 
strength discussed above is apparent also in its dynamics. 
In Fig. 0]the mean squared displacement is presented for 
different waiting times, as a function of t' = t — t w , for 
both quenches. As time proceeds, the system is more 
solid-like, and the displacements of particles is hindered. 
For comparison, (5r 2 (t')} for the equilibrium system at 
4> p = 0.42 is also presented. At 4> p — 0.80, even for short 
waiting times, the dynamics of the out-of-equilibrium 
systems is slower than the fluid one, although the mean 
number of bonds per particle is significantly smaller (see 
Figs. [T] and EJ. 

The comparison between both quenches shows that the 
system at <f) p = 0.80 slows down much faster than at 




t'/x 

FIG. 5: Incoherent density correlation function, $ a q (t), time 
rescaled to match the decay from the plateau. Different lines 
show different waiting times (same as Fig. , increasing from 
right to left. The thin green lines show for (f> p — 

(broken line) and 4>p — 0.42 (continuous line). 



4> p = 0.50. At large waiting times, the former shows a 
clear plateau in (6r 2 (t')) at ~ 6 x 10~ 3 a 2 , which extends 
for three decades in time. On the other hand, the state 
at (ftp — 0.50 only shows a hint of the plateau, which is 
located at longer distances, and the mean squared dis- 
placement reaches larger values. 

The self part of the intermediate scattering function, 
or incoherent density correlation function, shows simi- 
lar features as the mean squared displacement: a pro- 
gressive development of a plateau, clearer and higher for 
4>p = 0.80 than for p = 0.50. In Fig. [3 the correlation 
functions for qa — 9.9 (third peak in S(q)) are presented 
for the same waiting times as Fig. with the time axis 
rescaled in each case to data-collapse the final decay from 
the plateau. Very good scaling is observed for the initial 
decay from the plateau, whereas deviations are observed 
at long times. The scaling can be significantly improved if 
only the data with not-too-small waiting times are used, 
as shown below (Fig. I14fl . 

As expected, the plateau in the correlation function 
is clearer (flatter) in the case of 4> p = 0.80, whereas no 
real plateau is seen at <j> p = 0.50. When compared to 
4> p = (broken thin line), more than five decades in time 
of slowing down are observed. It is also noticeable that 
the shape of the relaxation curve changes with the at- 
traction strength; it is more stretched in the equilibrium 
system (<f> p = 0.42, continuous thin line) than in the aging 
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FIG. 6: Upper panel: Density correlation functions, coherent 
(thin lines, qa = 9.6) and incoherent (thick lines, qa = 15.6), 
for different states. From bottom to top: 4> p = 0.42 (fluid), 
0.50, 0.60, 0.80 and 1.00. Lower panel: Nonergodicity param- 
eters (coherent, open circles and incoherent, closed circles), 
and localization lengths vs. the distance to the transition 
point. The lines are square root fits to the data. 




qa 

FIG. 7: Upper panel: Nonergodicity parameter as a function 
of the wavevector for states progressively deeper in the glass. 
From top to bottom: critical, <j> p = 0.50, 0.60, 0.80 and 1.00. 
The inset presents the incoherent non-ergodicity parameter 
for the same states. Lower panel: Amplitude of the square- 
root fitting to f q , A q (black circles), and amplitude of the first 
order in the von Schweidler expansion of the a-decay, h q (red 
crosses) . 



states. We interpret these findings as arrested dynamics 
induced by particle bonding, indicated by the existence 
of a plateau in $ 9 for <fr p = 0.80, even with reversible 
bonds that are not stabilized artificially |2 1| . 

The evolution of the plateau can be observed in Fig. 
EJ where we study the (5 regime for <fi p = 0.50, 0.60, 
0.80 and 1.00. In this case, correlation functions were 
measured, for time intervals up to t' = 1000, for 100 
different waiting times in the range 10 4 < t w < 6 x 10 4 . 
(Six independent quenches were made for every state.) 
No evolution was observed, within the short times probed 
in Fig. El for either Q q (t) or the M.S.D., as expected from 
Figs. 01 and The height of the plateau in Q q (t) grows 
as the state is deeper in the glass, and the nonergodicity 
parameter is calculated from the value of the function in 
the (3 regime, i.e. between the two vertical dashed lines 
in Fig. Similarly, the localization lengths obtained 
from the plateau in the M.S.D. decreases as 4> p increases. 

Both the nonergodicity parameter (left scale in the 
lower panel, presented for the coherent, f q , and incoher- 
ent functions, /*) and the localization length, r/, show 
square root behaviours as a function of the distance to 

the transition: f q = f q + A q ^J <j) p — ^ . Here is the 

critical non-ergodicity parameter, A q is a g-dependent 



amplitude and <fr p is the transition point (with similar 
expressions holding for /* and r{). As an example, the 
lower panel of Fig. presents f q for qa — 9.6, for 
qa = 15.6 and r\ (right y-scale). Similar trends have been 
reported with molecular glass formers [3!|, and recently 
with m-Toluidine |40j , yielding a square root dependence 
on the temperature shift from the transition point. In- 
terestingly, MCT predicts such a dependence close to the 
transition on the glass side; but the square root regime 
observed here extends much further into the glass than 
does the asymptotic approximation. (The latter's range 
of validity can be found by numerical MCT calculations 
for our system |4l|.) 

The wavevector dependence of both the self and co- 
herent nonergodicity parameters is presented in Fig. [7| 
for the same states as Fig. Whereas f q shows oscilla- 
tions in phase with S(q), f q decays monotonically, both 
of which are in agreement with MCT expectations. For 
states deeper in the glass, the ^-distribution of f q or /* 
becomes wider, implying the decrease of the localization 
length. The critical nonergodicity parameter is obtained 
from the square-root fitting to the data, and agrees, to 
within our numerical accuracy, with results from equilib- 
rium simulations |15| . The critical (j> p for the transition 
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FIG. 8: Relaxation time normalized with Tt ra =8 as a function 
of the waiting time, for <f> p — 0.50 and 0.80, and the power-law 
fittings. The blue line is r ~ t w , i.e. simple aging. 



has been taken as 



t>p = 0.4265, in accordance with pre- 
vious results from equilibrium. The amplitude of the 



square-root in the fitting of f q vs. 



p 



also oscillates 



in (anti)phase with the structure factor (lower panel of 
Fig. UJ. Similar results were obtained with molecular 
glasses |3jj|40j, i.e., f q oscillates in phase with 5(g) and 
A q in antiphase; our results extend over the whole range 
of wavevectors. In ideal MCT, A q should be equal to 
the first order amplitude in the von Schweidler expan- 
sion of the structural a-decay, h q [3^|. The comparison 
of the two fitted quantities, shown in the lower panel of 
Fig. [7| is not satisfactory. This indicates that the fitted 
square-root dependence in Fig. is not the asymptotic 
MCT one; this concurs with the fact that the apparent 
square-root behavior extends over a much wider range of 
<f>p — cf>p than does the asymptotic MCT result. 

The relaxation time as a function of the waiting time 
is presented in Fig. [5] for both the <j> p — 0.50 and the 
4> p = 0.80 quench states. The ratio T tnl /T tw= s is given 
as a measure of the relaxation time; this is the necessary 
time factor in Fig. [3]to collapse all of the curves onto the 
correlator at t w = 8. Power laws are observed, r ~ 
with exponents bigger than one in both cases. Simple 
aging, for which r ~ t w , is also drawn in the figure; the 
fits show clearly that a > 1, so that the system evolves 
'faster' (in some sense) than its age. We might expect this 
behavior to cross over to linear aging eventually; mean- 
while it is reminiscent of "superaging" observed in spin 
glasses OH , although any mechanistic link between 
those systems and ours is rather unclear. 

The data collapse we observe under time rescaling of 
the correlation functions has also been seen in other 
structural glasses 2] , again accompanied by a power law 
relation between r and t w . In repulsion-driven glasses, 
this behavior sets in at t w = 0, but in colloidal gela- 
tion the power law dependence holds only at large waiting 
times It is preceded by a different behavior, 

crossing over either from another power-law pfi| , or from 
a more "explosive" regime (found also in experiments) 
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FIG. 9: Wavevector dependence of the relaxation time 
($£(t s ) = f q /2), for cj>p = 0.50 (black circles) and 0.80 (red 
squares) and t w — 4096, and for the equilibrium state at 
4> p = 0.42 for comparison (green diamonds). The dashed 
lines are power-law fittings with exponents (—2.48 ±0.02) for 
<p p = 0.50 and (2.25 ± 0.02) for <P P = 0.80, and (2.82 ± 0.03) 



for 



0.42 at high q and (1.92 ± 0.03) for low 



44, 45, 46]. In these cases, however, the kinetics of aging 
may be influenced by phase separation, which our study 
was specifically designed to avoid. 

The wavevector dependence of the time scale, defined 
via $ q (T q ) — f q /2, is presented in Fig. for both 
quench states at a waiting time t w = 4096 (chosen be- 
cause it shows a significant part of the structural relax- 
ation, allowing determination of T q ). In both states, the 
time scale follows a power law with q, with exponents 
greater than 2. In fluid states close to the glass transi- 
tion, according to MCT, this exponent is —1/6 where 
the von Schweidler exponent b describes stretching of 
the cv-decay. The results for the equilibrium state at 
4> p = 0.42 are presented for comparison, showing one ex- 
ponent (2.82 ±0.03) at high q (yielding b = 0.36, compat- 
ible with b — 0.37 obtained by direct von Schweidler fit- 
ting H3), and another exponent (1.92 ±0.03) for lower q, 
approaching the hydrodynamic limit. The exponents ob- 
tained here for our two quenches (around 2.4) are smaller 
than that of equilibrium states |25j , yielding higher val- 
ues of b, and hence less stretched decays of the density 
correlation function. This is indeed observed directly, in 
Fig. O Of the two quenches, the state at <j) p = 0.50 gives 
a higher exponent for the r q behavior, and its correlation 
function is more stretched indicating that the relation be- 
tween stretching and wavevector dependence of the time 
scale apparently holds in the glass. 

We study next the evolution of the bond correlation 
function with the waiting time. The bond correlation 
function, is defined as the fraction of bonds 

that have existed uninterruptedly since if = 0. (Re- 
call t' — t — t w .) This function is presented in Fig. 
for different waiting times, for both <j> p = 0.50 and 0.80. 
In both cases $s decays slower as the waiting time in- 
creases, indicating that the structure becomes stiffer, in 
agreement with the evolution of <I> g and of the structure 
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FIG. 10: Bond correlation function for the same waiting 
times as Fig. 0]from left to right and <j> p — 0.50 (upper panel) 
and 4> p = 0.80 (lower panel). Note the different y-scale in 
both panels. 



factors presented above. It should be noted also that the 
bonds are much more stable when the quench is made to 
4> p = 0.80, even at short waiting times, as expected from 
the preceding results. 

The bond correlation function decays in both states, 
what implies that bonds are never permanent, even at 
long waiting times (i.e. a plateau in &B(t') is not ob- 
served), although the fraction of bonds that remain un- 
broken is much larger for <f> p = 0.80 than <f> p — 0.50. In 
contrast, the mean squared displacement (Fig. an d 
density correlation function (Fig. [5J indeed showed ar- 
rested motion of the particles at long waiting times for 
(f) p = 0.80. These results from the bond correlation func- 
tion thus indicate that, although the relaxation of den- 
sity fluctuations and the diffusion of particles are both 
arrested, individual bonds are breaking (and new ones 
are being formed) within the cage or bonding network 
around each particle. 

This is further supported by analyzing the environ- 
ment, or cage, correlation function, $s(i'), which is de- 
fined as the fraction of bonds at time t' that existed at 
t' = (regardless of whether they broke in between). The 
function is presented for the same waiting times, in both 
quench states, in Fig. ^ As required, $£;(<') > $s(t') 
since bonds that break and reform are included in the 
former but not the latter. Note that for the shortest 
waiting time and 4> p = 0.50, less than 10% of the bonds 
have remained continuously unbroken, but, on average, 
more than half of a particle's bonds are still part of its 
environment. 

At long waiting times the environment correlation 
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FIG. 11: Same as Fig. 1101 for the cage correlation function 
(see text for the definition). Note that different y-scales are 
used in the two panels. 



function plateaus at a very high value, which extends to 
longer times as the waiting time grows. Thus, the net- 
work of particles comprising the system becomes stiffer 
and its relaxation takes longer. In the quench to <j) p — 
0.80 at the largest t w , has decayed only a few percent, 
implying that most of the particles have been unable to 
escape their network of bonds. 

These results for fully support the conclusions 
drawn from <I> g ; in fact, <&e could almost be viewed as a 
density correlation function at a wavevector correspond- 
ing to the size of the neighborhood of the particles. In- 
terestingly, some bonds can break even at long waiting 
times (as shown in Fig. 110(1 , although the local structure 
is frozen. Our results show that the dynamical arrest is 
a collective phenomenon that can be induced by bonds 
with finite energy, and thus finite life-time, in contradic- 
tion with previous works where it was concluded that the 
attractive glass was unstable due to the finite bonding 
energy pif . 

To further exploit the relationship between structural 
arrest and bond formation, we time-rescale the bond cor- 
relation function and density correlation function $ g 
using the waiting time, t w ; see Fig. ^] for <f> p = 0.50. 
After the rescaling, neither $ g or <E>b collapse, but the 
density correlation functions cross at x — t/t w « 10 , 
whereas no crossing is observed for $s(x). Therefore, 
at constant x beyond the crossing point, i.e. x > 10 , 
$*(x) is larger for longer waiting times, whereas $.b(x) 
is smaller. In other words, the time scale for the slow- 
ing down of structural relaxation grows faster than the 
'bond life time'. (Note that a unique bond life time can- 
not be defined, however, because the decays of the bond 
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FIG. 12: Density correlation lunction (upper panel) and bond 
correlation function (lower panel) time rescaled with waiting 
time, t w , for <f> p = 0.50. Different colors correspond to differ- 
ent waiting times, as in Fig. 4. 
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FIG. 13: Upper panel: Structure factor of the system without 
repulsive barrier for different waiting times, as labeled. Lower 
panel: Comparison of the structure factors with the same 
waiting time, t w = 4096, for the states <j> p = 0.50, (f> p = 0.80, 
and <j) p = 0.80 without the repulsive barrier, as labeled. 



correlation functions do not collapse for different waiting 
times.) This decoupling between structural relaxation 
and bond formation again points to the collective origin 
of the attraction-driven glass transition. The implica- 
tions of Fig. 1121 differ from previous results for equi- 
librium systems where the dynamics (as measured by 
the mean squared displacement) and the bond correla- 
tion function were scaled using thermodynamic proper- 
ties (the second virial coefficient, B2) 26|. In our case, 
no quantity can be used collapse by time-rescaling both 
$ ? and <&b, because they evolve differently. Currently 
we are unable to reconcile our findings with those of Ref. 
pjfjj. if the latter are to apply at low temperatures (high 
interaction strength). We can only suspect that there is 
some fundamental difference in either the physical speci- 
fication of the system or the computational methodology. 
However, we have not been able to identify such a differ- 
ence. 

It is also worth noticing in Fig. ^| that the rescaled 
curves, both for <I> 9 and <&b, approach each other more 
closely at long waiting times. It is tempting to expect col- 
lapse at still longer waiting times, implying a crossover 
from the "superaging" regime reported above to a simple 
aging one, thus leaving t w as the only physical time scale 
in the system. Our present simulations, however, only 
hint at this possibility, which we can neither prove nor 
rule out. Similar results are obtained for <p p — 0.80 (data 



not shown), both in the decoupling of the time scale for 
structural relaxation and bond life time, and in the sug- 
gested simple aging regime for long waiting times. 



C. Effect of the repulsive barrier 

In order to test the effect of the repulsive barrier in 
the interaction potential, introduced to prevent phase 
separation but possibly inducing microphase separa- 
tion, we have performed simulations without the barrier. 
The system, thus, is closer to the experimental colloid- 
polymer mixtures where gels have been deeply studied 
pi l34l l3(| . In previous simulation works, it was shown 
that the liquid-gas separation is arrested by the attrac- 
tion driven glass transition, in agreement with experi- 
ments [23, EI H2- We show below that indeed this is 
the case in our system and we analyze the effect of the 
repulsive barrier on the dynamics of the system. 

In Fig. ^|we show the evolution of the structure fol- 
lowing the quench. The structure factor develops a low-q 
peak which moves to lower wavevector and grows, but ap- 
parently saturates during the simulation time - for com- 
parison, in the system with <f> p = 0.35 below the glass 
point cj)p = 0.4265 the low-q peak moves out below the 
accessible wavevectors, for t w > 128. The structure fac- 
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FIG. 14: Master decay of the incoherent density correla- 
tion function for cj> p — 0.80 without repulsive barrier (thick 
black lines) and with it (thin red lines) and <f> p = 0.50 with 
it (thin green lines). The master decays have been obtained 
from systems with t w = 256, 512, 1024, 2048, 4096, 8192, 
16384, from right to left in upper part. The inset shows the 
non-ergodicity parameters for the same three states, and the 
critical one (from Ref. @), as a function of the wavevector. 



tors of states with and without the repulsive barrier are 
compared for the same waiting time in the lower panel of 
Fig. The peak is placed at lower wavevectors and is 
higher for the system without barrier. These results indi- 
cate that the phase separation is prevented by the glass 
transition, resulting in a heterogeneous system with voids 
and dense regions. The effect of the repulsive barrier in 
the potential is then to reduce the size and concentration 
of voids. 

An analysis similar to that of Sec. III.B can be con- 
ducted to study the dynamics of the arrested system. 
The mean squared displacement and density correlation 
function show a slowing down of the dynamics as the 
waiting time increases, in a two-step process, similar to 
glass aging. The density correlation functions can be 
rescaled to overlap the decay from the plateau as shown 
in Fig. ^]for <f> p = 0.80 without barrier, where only wait- 
ing times above t w = 256 are included. For comparison, 
the master function for <j> p = 0.80 and <j) p = 0.50 with the 
repulsive barrier are also included. Note that the scaling 
of the curves is much better than in Fig. [5] since only 
data with t w > 256 is used. The two master functions 
for <fip = 0.80 are very similar over the whole decay from 
the plateau, indicating that the long time relaxation of 
the system is only slightly affected by the repulsive bar- 
rier. The time scale of this relaxation as a function of the 
waiting time (not shown) also obeys a power law, with 
exponent larger than one: t ~ t^ 20 , very close to the 
exponent found in the system with barrier, 1.16. 

The non-ergodicity parameters for the incoherent cor- 
relation function are presented in the inset of Fig. 1 141 for 
the same states. The agreement between the two /| at 
4> p = 0.80 indicates that the localization length is similar 
in both cases, and that the origin of the transition is the 
same: bond formation. Fig. 1151 presents the evolution 




FIG. 15: Bond correlation function for the system with <j> p = 
0.80 without the repulsive barrier, for the same waiting times 
as Fig. [I] (upper panel) . The lower panel shows a comparison 
of the bond correlation function for systems at cj> p = 0.80 
without the barrier (thick lines) and with the barrier (dashed 
lines) for t w = 64, 1024 ad 16384 (from bottom to top). 



of the bond correlation function for the state P = 0.80 
without the repulsive barrier (upper panel). This func- 
tion decays even at very long waiting times, although the 
mean bond-life time grows continuously. The compari- 
son of the bond correlation functions for the system with 
and without the barrier (lower panel), shows completely 
similar trends, and that the bonds for the system without 
the barrier are slightly stronger. 

The results presented in this subsection indicate that 
the repulsive barrier plays a minor role in the behaviour 
of the system beyond the glass transition, although it has 
a tremendous effect below it suppressing the liquid-gas 
separation. The competition between glassy physics and 
phase transition results in a system with arrested phase 
separation, whose properties are largely similar to those 
of the system without the phase transition. This study 
also rules out the possibility that microphase separation, 
if it occurs in our system, is primarily responsible for the 
results presented above. 



IV. CONCLUSIONS 

In conclusion, we have shown that aging of attraction 
driven glasses is broadly analogous to the aging of re- 
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pulsion driven glasses, albeit with some specificities aris- 
ing from the different driving mechanism in either case: 
bond formation for the former, steric hindrance for the 
latter. Specific differences appear, for instance, in the 
values of localization length and in f q . In attraction 
driven systems, as in repulsive ones, the mean squared 
displacement and density correlation functions display a 
plateau followed by a decay, which moves to longer times 
as the waiting time grows. $q(t') can be time-scaled 
to collapse during the decay from the plateau for differ- 
ent waiting times, and the decay time follows a power- 
law dependence on the waiting time, with an exponent 
slightly bigger than unity. The non-ergodicity parameter 
increases as the state is deeper in the glass, following an 
apparent square-root behavior, and oscillates in phase 
with the structure factor, decaying at larger wavevec- 
tors with increasing <p p , implying a shortening localiza- 
tion length. The bond and environment correlation func- 
tions show that as the aging proceeds and more bonds are 
established, a stiffer network of particles forms whose re- 
laxation time increases; but that the bonds themselves 
are not permanent. However, whereas the time scale for 
structural relaxation grows faster than the waiting time, 
the bond life time increases slower than t w . This decou- 
pling of structural relaxation from single bond breaking 
illustrates the collective character of the dynamical ar- 
rest. Yet, from our simulations, it cannot be ruled out 
that for large enough t w both time scales approach simple 
aging behavior, scaling directly with t w . The effect of the 
repulsive barrier (introduced to suppress liquid-gas sep- 
aration) on the dynamics is quite weak beyond the glass 
transition; this suggests that the transition is primarily 
induced by bonding, and that aging in attractive glasses 
is dominated by the evolution of the bonds rather than 
by any slow evolution of the phase separation pattern. 

Important differences can be observed between the 
quenches to <f> p = 0.50 and <j> p = 0.80. The former 
is slightly above the transition line, and does not show 
clear plateaus in <!> g , the mean squared displacement or 
$£. The stretching of the correlation functions is also 
different (more stretched closer to the transition point). 
Nevertheless, scaling of the density correlation functions 
was fulfilled, and the exponent of the relaxation time as 
a function of the waiting time is similar for both states. 



From the structural point of view, the glass with lower at- 
traction strength allows for more restructuring: the holes 
are bigger, as noticed by a higher low-g peak in the struc- 
ture factor, and more bonds per particle are formed. In- 
terestingly, our molecular dynamics results strongly in- 
dicate that bonds of finite energy, which are therefore 
reversible on the time scale of the simulation, can induce 
dynamical arrest at high enough interaction strength. 
These conclusions are different from those drawn from 
previous works [24| although the results themselves are 
not very different from ours. The decoupling of structural 
relaxation from single bond breaking represents another 
aspect where our results differ from previous findings |26j | , 
if those can be extrapolated to low temperatures. 

The softness of the glass at low attraction strength 
(close to the boundary, <p p — 0.50) may arise from a non- 
vanishing population of highly mobile particles detected 
in the fluid phase close to the transition point 0, • 
This population coexists with a rigid network of particles, 
and reduces as the attraction strength is increased, being 
completely absent at (f> p = 0.80 28]. The origin of the 
softness can thus perhaps be traced back to the strong 
dynamical heterogeneities in the fluid. 

Finally, in addition to the results presented in this pa- 
per, we have found recently that many dynamical quan- 
tities for the fluid states close to the glass line can be 
accounted for within MCT, not only in its universal qual- 
itative predictions, but also quantitatively fl5L l4l| . Fur- 
thermore, the dynamic heterogeneities of the glass de- 
crease as the state point moves deeper into the attractive 
glas s | 2 8l , in agreement with findings for repulsive glasses 
\VA Il8j . All of these results support the argument that 
attraction driven glasses share many universal properties 
with repulsion driven glasses, albeit with some differences 
that are attributable to the different driving mechanism. 
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